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Abstract 

Background: Endogenous small (sm) RNAs (primarily si- and miRNAs) are important frans/c/s-acting regulators 
involved in diverse cellular functions. In plants, the RNA-dependent RNA polymerases (RDRs) are essential for smRNA 
biogenesis. It has been established that RDR2 is involved in the 24 nt siRNA-dependent RNA-directed DNA 
methylation (RdDM) pathway. Recent studies have suggested that RDR1 is involved in a second RdDM pathway that 
relies mostly on 21 nt smRNAs and functions to silence a subset of genomic loci that are usually refractory to the 
normal RdDM pathway in Arabidopsis. Whether and to what extent the homologs of RDR1 may have similar 
functions in other plants remained unknown. 

Results: We characterized a loss-of-function mutant {Osrdrl) of the OsRDRl gene in rice {Oryza sativa L.) derived 
from a retrotransposon Tos17 insertion. Microarray analysis identified 1,175 differentially expressed genes (5.2% of all 
expressed genes in the shoot-tip tissue of rice) between Osrdrl and WT, of which 896 and 279 genes were up- and 
down-regulated, respectively, in Osrdrl. smRNA sequencing revealed regional alterations in smRNA clusters across 
the rice genome. Some of the regions with altered smRNA clusters were associated with changes in DNA 
methylation. In addition, altered expression of several miRNAs was detected in Osrdrl, and at least some of which 
were associated with altered expression of predicted miRNA target genes. Despite these changes, no phenotypic 
difference was identified in Osrdrl relative to WT under normal condition; however, ephemeral phenotypic 
fluctuations occurred under some abiotic stress conditions. 

Conclusions: Our results showed that OsRDRl plays a role in regulating a substantial number of endogenous 
genes with diverse functions in rice through smRNA-mediated pathways involving DNA methylation, and which 
participates in abiotic stress response. 
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Background 

RNA silencing is an evolutionally conserved gene regula- 
tion mechanism in eukaryotes mediated by 20-25 nt non- 
coding small (sm)RNAs. These smRNAs are processed 
from double-stranded (ds) or hairpin RNA molecules by 
Dicer-like (DCL) proteins, and guide RNA-induced silen- 
cing complexes to cognate single-stranded RNAs based on 
sequence complementarity, and result in degradation of 
the targeted RNAs [1-3]. In plants, there are several differ- 
ent classes of smRNAs, including 20-24 nt micro RNAs 
(miRNAs) processed by DCL1, 21-22 nt small interfering 
RNAs (siRNAs) by DCL4 and DCL2, and the 24 nt 
heterochromatin- associated siRNAs by DCL3. miRNAs 
play an important role in plant development by directing 
posttranscriptional gene silencing (PTGS) of regulatory 
genes such as those encoding transcription factors. 
Similarly, 21-22 nt siRNAs guide the degradation of viral 
RNAs as well as some endogenous mRNAs and are im- 
portant for plant defense against viruses and for some 
aspects of plant development [4-6]. Unlike these PTGS- 
associated smRNAs, the 24 nt siRNAs are associated 
with RNA-directed DNA methylation (RdDM), a plant- 
specific de novo DNA methylation pathway required for 
transcriptional silencing of transposable elements and 
other DNA repeats to maintain genome stability [7-10]. 

The biogenesis of siRNAs in plants requires the activ- 
ity of RNA-dependent RNA polymerase (RDR), which 
converts single-stranded RNAs to dsRNA precursors of 
siRNAs. The dicot model plant Arabidopsis thaliana 
has six RDR genes, i.e., RDR1, RDR2, RDR3a, RDR3b, 
RDR3c and RDR6 [11], of which three RDRs (RDR1, 
RDR2 and RDR6) are shown to play roles in the RNA 
silencing pathways. RDR2 is required for 24 nt siRNA 
biogenesis and therefore involved in the canonical 
RdDM pathway [7-9]. RDR6 is involved in the produc- 
tion of the endogenous 21 nt trans-acting siRNAs and 
also essential for sense transgene-induced PTGS 
[12,13]. Both RDR6 and RDR2 are also involved in viral 
siRNA accumulation in infected Arabidopsis plants 
[14-16]. The function of RDR1 in RNA silencing is less 
understood, but recent studies have shown that it is 
involved in siRNA biogenesis from a subset of RNA vi- 
ruses [17-19]. Furthermore, rdrl mutant of Arabidopsis 
showed loss of DNA methylation in a subset of genomic 
loci in comparison to wild-type Arabidopsis plants [20,21], 
suggesting that RDR1 plays a role in the recently identified 
non-canonical, 21 nt siRNA-directed RdDM pathway 
[20,21]. However, the function of RDR1 in gene regulation 
from a genome-wide perspective has not been investigated 
in any plant. 

In contrast to Arabidopsis that has six RDR genes, the 
RDR family of rice (Oryza sativa L.), a model plant for 
monocots, contains only three members, namely OsRDRl, 
OsRDR2 and OsRDR6 [11,22,23]. A previous study showed 



that OsRDRl has a similar function to its counterparts 
in Arabidopsis and tobacco (Nicotiana tabacum) in PTGS- 
based silencing of certain RNA viruses, such as Bromovirus 
[6,19,24]. To investigate if OsRDRl plays a role in regula- 
tion of endogenous genes in rice, we characterized a loss- 
of-function mutant of OsRDRl derived from a disruptive 
LTR retro transposon (ToslT) insertion into the 2nd exon 
of the gene. We investigated genome-wide changes in 
gene expression and smRNA profiles, localized changes 
in DNA methylation, and phenotypes under normal and 
several abiotic stress conditions in this rice rdrl mutant. 

Results 

Characterization of the rice RDR1 mutant (Osrdrl) 

We obtained a LTR retro transposon Tosl7 [25] insertion 
line for OsRDRl (accession number H0643) from the 
Tosl7 insertion mutant library of rice cv. Hitomebore 
(www.cns.fr/spip/Oryza-sativa-retrotransposon-Tosl7.html). 
Molecular characterization identified H0643 as hetero- 
zygous for a Tosl7 insertion into the second exon of 
OsRDRl (Figure la). We obtained the homozygous mutant 
(OsRDRl-/- or Osrdrl) and its sibling wild type (WT) 
plants by selfing of the heterozygous plant (OsRDRl+/-) 
for five successive generations. In each generation, the 
three kinds of genotypes, WT, heterozygote and homo- 
zygous mutant, were selected based on locus-specific 
PCR amplifications (Figure la). Both the heterozygous 
and homozygous plants for OsRDRl showed no discernibly 
altered phenotypes in the entire growth and developmental 
period over multiple generations under normal field condi- 
tions (Figure lb). Semi-quantitative and real-time quantita- 
tive (q)RT-PCR analyses confirmed that the homozygous 
OsRDRl mutant (Osrdrl) had a complete loss of OsRDRl 
expression in shoot-tip tissue wherein the gene was highly 
expressed in WT plants (Figure lc). This indicated that 
the exonic insertion of Tosl7 knocked out the expression 
of OsRDRl, and hence, abolished its function. 

Genome-wide changes in gene expression in Osrdrl 

We profiled the transcriptome of shoot-tip tissues between 
Osrdrl and its sibling WT plants using the Affymetrix Gen- 
eChip Rice® Genome Array (The Affymetrix, Inc. Santa 
Clara, CA, USA). After normalization of the microarray 
data, we detected 57,381 expressed genes in the shoot- tip 
tissue of rice. The expression levels of 22,419 genes were 
conserved between Osrdrl and WT, but 896 and 279 
genes showed significant up- and down-regulation in 
Osrdrl, respectively (Figure 2a). A Gene Ontology (GO) 
category analysis of these 1,175 differentially expressed 
genes showed that they were enriched in a variety of GO 
categories (Figure 2b). However, these 1,175 differentially 
expressed genes were found to distribute non-randomly 
across the 12 rice chromosomes (P = 2.2E-16, based on 
Chi square test). For example, chromosomes 1, 2 and 3 
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Figure 1 OsRDRl gene expression is abolished in the 7bs77 insertion mutant Osrdrh (a) Structure of the OsRDRl locus with the Tos17 
insertion into the 2nd exon (vertical arrows). Heterozygous/homozygous individuals were selected based on PCR which primers were indicated 
by the purple arrows while primers for OsRDRl expression analysis were indicated by black arrows, (b) Germinated seedling, paddy-field-grown 
plant and kernel phenotypes of wide-type (WT) and Osrdrl. (c) Relative to WT, OsRDRl expression was silenced in shoots of Osrdrl, as evidenced 
by both qRT-PCR (top) and semi-quantitative RT-PCR (bottom) amplifications with gene-specific primers downstream of the Tosl7 insertion 
(horizontal arrows). Genomic DNAs were used as positive controls. 



contained significantly more distributions than the rest 
chromosomes (Figure 2c). It is also clear from the data 
that within a given chromosome, the distribution is also 
nonrandom, for example, the distributions are almost ex- 
clusively confined to the long arms of chromosomes 4 and 
8 and 9 relative to their respective short-arms (Figure 2c). 

The highly reproducible microarray profiles among 
three biological replicates for both Osrdrl and its WT 
sibling plants testified reliability of the data and their ana- 
lysis. All microarray data have been submitted to the GEO 
repository under the accession number of GSE58007. To 
further verify the quality of the microarray data and 
analysis, we analyzed 18 genes representing both up- 
and down-regulation in Osrdrl vs. WT, as well as equal 
expression between the two lines using qRT-PCR assay 
on the same cDNAs as used for microarray. The qRT- 
PCR results were highly consistent with the microarray 
data for almost all the 18 tested genes in levels or at 
least in trends of expression changes (Figure 2d), con- 
firming reliability of the microarray analysis. 

Alteration in smRNA clusters in Osrdrl 

Previous studies have established that RDR1 function is 
required for biogenesis and/or amplification of some 
types of RNA virus-related smRNA accumulation in 
Arabidopsis [26] and tobacco {Nicotiana tabacum) [27]. 
These findings promoted us to test whether loss of func- 
tion of OsRDRl may have a general impact on "normal" 
smRNA abundance in rice, and we investigated this issue 



by high-throughput smRNA sequencing. Comparison of 
the 10,398,592 clean smRNA reads from Osrdrl with 
the 9,339,435 reads from its sibling WT {see Methods) 
revealed highly similar profiles in both size distributions 
and sequence categories of the smRNAs between Osrdrl 
and WT (Figure 3a, b), suggesting that the overall smRNA 
abundance was not generally affected by the loss of func- 
tion of OsRDRl. 

Genome-wide overall similarity in smRNA abundance 
does not necessitates absence of smRNA fluctuations 
in localized smRNA clusters, because up- and down- 
regulated smRNA accumulation can be masked by re- 
ciprocal compensation. We thus investigated localized 
smRNA accumulation between Osrdrl and its sibling 
WT by mapping the cleaned smRNA reads to 100 bp 
sliding windows (being reflected as smRNA clusters) 
across each of the 12 rice chromosomes, normalized the 
reads to Reads per Million (RPM), and then compared the 
RPM smRNA clusters between Osrdrl and WT. Using 
4-fold difference as a cut-off threshold, we identified 
many smRNA clusters with differential abundance be- 
tween Osrdrl and its sibling WT, which were uniformly 
distributed across the entire length of each chromosome 
(Figure 3c). Next, we extracted the differentially expressed 
smRNAs between Osrdrl and WT (also based on a cut-off 
threshold of 4-fold difference) in the size range of 
20-24 nt, which should parsimoniously contain all siRNAs, 
and mapped them to the same 100 bp windows across each 
chromosome. We found that this subset of differential 
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Figure 2 Effects of null mutation of OsRDRI on global gene expression in rice, (a) A summary of microarray analysis showing the total 
number of genes detected and the number of differentially expressed genes in WT and Osrdrl (b) Gene Ontology (GO) category analysis of the 
1,175 differentially expressed genes between WT and Osrdrl . The y-axis is the percentage of genes mapped by the GO category terms: the 
percentages were calculated by the number of Osrdrl vs. WT differentially expressed genes divided by the total number of genes mapped to the 
particular GO category. The x-axis is the GO category terms which were ordered by their relative abundance (total number of expressed genes in 
the shoot-tip tissue of rice are 45,078). The blue bars denote percentages for each category of all the annotated genes, and the red bars denote 
percentages of the GO categories of all the expressed genes, (c) Distribution of the 1,175 Osrdrl vs. WT differentially expressed genes across each 
of the 12 the rice chromosomes (horizontal lines). The red columns above and the green columns below the chromosomes represent up- and 
down-regulated genes in Osrdrl vs. WT, respectively. The y-axis indicates fold changes in gene expression between Osrdrl and WT. The vertical 
blue bars denote for centromeric regions, (d) Validation of the microarray results by qRT-PCR, where the blue and red columns represent WT and 
Osrdrl, respectively. Numbers 1 and 2 below each gene represents data of microarray and qRT-PCR results, respectively. Statistical significance at 
that P <0.05 and P <0.01 levels is marked by one or two asterisks. 



smRNA clusters also distribute on both arms of each 
chromosome (Additional file 1), although due to their 
smaller numbers, we cannot rule out the possibility that 
the distribution might show "hot spots" within a given 
chromosome. Taken together, the smRNA sequencing 
data suggested that loss-of-function mutation in OsRDRI 



caused extensive alterations in smRNA clusters across 
each chromosome and throughout the genome, but it did 
not result in marked fluctuations of overall smRNA pro- 
files, probably due to more or less equally increased and 
decreased abundance of the smRNA clusters which offset 
each other. 
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Figure 3 Effects of null mutation of OsRDRI on genome-wide profiles of smRNAs in rice, (a) Size distribution of smRNA from Osrdrl and 
WT; (b) Categories of smRNAs from Osrdrl and WT; the red and blue columns represent Osrdrl and WT, respectively, (c) Difference of smRNA 
clusters (RPM) within 100 bp sliding windows across each of the 12 rice chromosomes between Osrdrl and WT, where x axis is the length of 
chromosome and y axis is the value of different RPMs (log value, base 2). The vertical blue lines denote centromeric regions in each chromosome. 



Altered expression of miRNAs and their target genes in 
Osrdrl 

Given the diverse important roles played by miRNAs, 
we investigated if their accumulation might be affected 
in the Osrdrl mutant Previous computational and cloning 
studies have identified ca. 300 miRNAs from 86 miRNA 
families in rice [28,29]. Based on this information, we first 
analyzed the abundance of known rice miRNAs (Osa- 
miRNAs) (listed in miRBasel4.0) in Osrdrl and WT. 
This analysis (Figure 4a, b; Additional file 2) indicated 
that: (1) majority (90.7%) of the known Osa-miRNAs 
were expressed equally or nearly so between Osrdrl and 
WT; (2) some miRNAs (5%) showed > 2 fold increased 



expression in Osrdrl relative to WT, with the highest 
expression ratio reaching 9.0:1 (£-test, P < 0.05); (3) 
some miRNAs (4.3%) showed >2 fold decrease in expres- 
sion in Osrdrl relative to WT, with the lowest expression 
ratio of 1:6.1 observed for miR167j between mutant and 
WT (£-test, P < 0.05); (4) Osa-miR395p and Osa-miR395s, 
being from the same miRNA family, showed changes in 
expression to opposite directions, with an expression ratio 
of 1:5.6 for osa-miR395p but 4.2:1 for osa-miR395s in 
mutant vs. WT (£-test, P < 0.05). To verify the expression 
differences based on the smRNA sequencing data, we per- 
formed semi-nested qRT-PCR analysis of four miRNAs in 
mutant and WT. The qRT-PCR results were found 
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Figure 4 Effects of null mutation of OsRDRl on abundances of miRNAs in rice, (a) Expression comparison of known miRNAs between 
Osrdrl and WT. (b) Three patterns of expression changes of known miRNAs. The red and blue spots represent up and down-regulated in Osrdrl, 
respectively; the green spots represent no variation between Osrdrl and WT. (c) qRT-PCR to verify variation of miRNA accumulation between 
Osrdrl and WT. Statistical significance is marked by asterisks. 

V J 



consistent with the smRNA sequencing data (Figure 4c), 
confirming the changes in miRNA expression between 
Osrdrl and WT. 

In addition to the known miRNAs, we identified a 
total of 10 putative novel osa-miRNAs from the smRNA 
data of the mutant and WT plants based on prediction 
of pre-miRNA-like stem-loop structures in sequences 
surrounding the smRNA sequences in the rice genome 
(Additional file 3). Three of these novel miRNAs (Osr- 
miRNA-N5.1, -N5.2 and -N5.3) had an identical mature 
sequence but corresponded to three independent genomic 
loci, thereby forming a novel miRNA family. For the puta- 
tive miRNA Osr-miRNA-N7, smRNA reads were detected 
from both the 5' (5p) and the 3' (3p) half of the predicted 
stem-loop structure, but the 5p smRNA showed a higher 
abundance than the 3p smRNA (Additional file 3), indicat- 
ing that the 5p smRNA is the guide strand (miRNA) 



whereas the 3p smRNA is the passenger strand (miRNA*) 
[30]. Like the known miRNAs, these novel miRNAs also 
showed expression variation between Osrdrl and WT, 
with 3 showing expression only in Osrdrl, and 5 showing 
expression only in WT plants (Additional file 3). Taken 
together, our results suggest that OsRDRl was likely 
involved in miRNA accumulation in rice. 

To investigate if the altered miRNA accumulation in 
Osrdrl relative to WT was associated with changes in 
miRNA target gene expression, we compared the miRNA 
expression profiles (abundance) derived from the smRNA 
sequencing data with the target gene expression levels 
based on the microarray data. We did not find a general- 
ized relationship between the miRNA abundance and tar- 
get gene expression levels (Additional file 4a). Instead, 
four types of relationships were recognized for a subset of 
miRNAs and their predicted targets (Additional file 4b), 
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which included: (1) reduced miRNA abundance was 
correlated with up-regulated expression of target genes 
in Osrdrl relative to WT (Additional file 4b-i); (2) in- 
creased miRNA abundance was correlated with down- 
regulated expression of target genes in Osrdrl relative 
to WT (Additional file 4b-ii); (3) both miRNAs and 
their target genes were up-regulated in Osrdrl relative 
to WT (Additional file 4b-iii); (4) both miRNAs and 
their target genes were down-regulated in Osrdrl rela- 
tive to WT (Additional file 4b-iv). The first two types of 
relationships supported a role of miRNAs in down- 
regulating expression of their predicted target genes. 
The last two types of relationships could be a result of 
concordant transcriptional regulation of the miRNAs 
and their target genes caused by another more upstream 
regulator(s) whose expression or activity was modified 
due to loss of function of OsRDRl. All small RNA data 
have been submitted to GenBank under the accession 
numbers of SRP042238. 

Locus-specific alteration of DNA methylation in Osrdrl 

As the Arabidopsis RDR1 has been shown to play a role 
in the non-canonical, NERD-dependent RdDM pathway 
[20,21], we were interested to know if OsRDRl might 
have a similar function in rice. We therefore examined 
cytosine methylation and gene expression levels of 10 se- 
lected genomic loci in Osrdrl and its sibling WT using 
bisulphite sequencing and qRT-PCR analysis. These loci 
overlapped with two transposable elements (TEs) and 
three protein-coding genes, which were chosen as repre- 
sentatives because they all showed alteration in smRNA 
clusters in Osrdrl relative to WT (Additional file 5). The 
bisulphite-sequenced regions for the two TEs (retro- 
transposon Tosl7 and DNA transposon Pong) included: 
(1) portions of the 5'- and 3'-LTRs together with their 
immediate flanking regions of two Tosl 7 copies (located 
on chromosomes 10 and 7, respectively) (Additional 
file 5a, b); (2) the 5' termini along with their immediate 
flanking regions of two Pong copies (located on chro- 
mosomes 2 and 9, respectively) (Additional file 5c, d), 
and; (3) a body-region of the transposase-encoding ORF 
of Pong (Additional file 5e) that is shared by all conserved 
copies of the element. The bisulphite-sequenced regions 
of the three protein-coding genes are all within their 
5'-upstream regions (Additional file 5f, g, h). 

The bisulphite sequencing results showed that: (1) of 
the five Tosl 7 regions analyzed, only the 5' LTR region 
for the Tosl 7 copy located on chromosome 7 showed 
marked decrease (by ca. 30%) in CG and CHG methy- 
lation but not in CHH methylation in Osrdrl relative 
to WT (Figure 5a and Additional file 5b); (2) for the 
three Pong regions analyzed, only the 5'region of the 
copy located on chromosome 2 showed clear methyla- 
tion changes: decrease in CG methylation by 20% and 



increase in both CHG and CHH methylation by approxi- 
mately 30% and 50%, respectively, in Osrdrl relative to 
WT (Figure 5a and Additional file 5c); of the three genie 
loci analyzed, only one (Os06g03 16000) showed increase 
in CG methylation by ca. 50% in Osrdrl relative to 
WT, while methylation of the other two regions were 
unchanged (Figure 5a and Additional file 5f). 

The two TEs showed significant up-regulation in Osrdrl 
relative to WT (Figure 5b), consistent with a decrease in 
methylation at the 5' regions of one copy of each TEs in 
Osrdrl (Figure 5a). Notably, all the three genes analyzed 
did not show the expected relationship between DNA 
methylation state of their 5'-regulatory regions and expres- 
sion levels. Specifically, one gene (Os06g03 16000) that 
showed increase in CG methylation in Osrdrl was up- 
regulated in expression (Figure 5b); the remaining two 
genes showed down-regulation in Osrdrl relative to WT 
despite the lack of methylation changes in the bisulphite- 
sequenced regions (Figure 5b). 

We next investigated possible relationships between 
smRNA accumulation and DNA methylation. We found 
that almost all of the altered CHH methylation was associ- 
ated with changes in smRNA clusters. For example, the in- 
creased CHH methylation of the Pong copy located on 
chromosome 2 was associated with a moderate increase in 
smRNA accumulation, whereas the slight decrease in 
CHH methylation in the flanking region and the increase 
in CHH methylation in the gene body region of the Pong 
copy located on chromosome 9 were associated with 
moderate decrease and increase in smRNA accumula- 
tions, respectively (Additional file 5c). These positive 
correlations of smRNA accumulation and CHH methy- 
lation suggests that OsRDRl plays a role in the de novo 
CHH methylation in a subset of genomic loci in rice, 
probably by affecting the production/accumulation of 
smRNAs required for RdDM, as shown in Arabidopsis 
[20,21]. The locus-specificity of methylation changes or 
the two analyzed TEs indicated that their methylation 
patterns were determined by either or both the flanking 
sequences and the local chromatin environment, an 
issue which warrants further investigations. 

Phenotypes in Osrdrl under normal and abiotic stress 
conditions 

It is known that various stress conditions may produce 
protracted effects on genome stability, leading to trans - 
generational changes in genome structure, which are 
proposed to have been initiated by epigenetic mecha- 
nisms [31-37]. We were therefore interested to know if 
OsRDRl may play a role in stress response in rice. We 
quantified phenotypes between Osrdrl and WT plants 
under normal and several short-term abiotic stress con- 
ditions {see Materials and Methods), which included 
treatments with salt, heavy metals Cu 2+ and Hg 2+ , and 
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Figure 5 Effects of null mutation of OsRDRI on locus-specific alterations of DNA methylation in rice, (a) Alteration in DNA methylation 
between Osrdrl and WT in the three cytosine sequence contexts, CG, CHG and CHH, based on bisulphite sequencing of 10 genomic loci from 
two transposable elements (TEs), Tosl 7 (four regions) and Pong (three regions) and three genes (one region each), (b) Expression differences of 
the TEs and genes between Osrdrl and WT based on qRT-PCR analysis. 



overdose nitric oxide (NO). The results showed that no 
phenotypic difference was found between Osrdrl and WT 
under normal condition, but significant ephemeral pheno- 
typic differences between the two genotypes emerged in 
some of the different stress conditions (Figure 6). Specific- 
ally, (1) seedlings of Osrdrl were more sensitive than WT 
to salt and overdose NO treatments, as being reflected by 
reduced plant height, root length and biomass at the seed- 
ling stage, with the difference in plant height and root 
length being persisted to the heading stage after removal 
of the stresses; (2) Seedlings of Osrdrl showed increased 
tolerance to heavy metal Cu 2+ /Hg 2+ as indicated by in- 
creased root length, but upon removal of the stresses the 
differences were gradually attenuated and completely dis- 
appeared at the heading stage; (3) When both unstressed 
and the transiently-stressed plants of the mutant and WT 
were grown to maturity, no difference in plant height, 
tiller number, panicle and kernel traits was observed be- 
tween the two genotypes. Collectively, our results suggest 
that OsRDRI has a potential function in stress response in 
rice, but the effects are contingent with presence of 
stresses without exerting protracted influence when the 
stresses are removed. 



Discussion 

RNA silencing pathways have been well characterized in 
the dicot model plant Arabidopsis but remain poorly stud- 
ied in other plants like monocots to which many major 
crops belong. Utilizing a retransposon Tosl 7 insertion 
mutant of OsRDRI in rice, we performed genome-wide 
analysis to unveil the function of OsRDRI, a component 
recently shown in Arabidopsis to be involved in non- 
canonical, 21 nt siRNA-directed RdDM pathway [20,21], 
on expression of endogenous genes. By deep sequencing 
of smRNAs and microarray analysis of gene expression in 
the Osrdrl mutant and its sibling WT, we showed that the 
expression of > 1,000 genes were significantly changed in 
Osrdrl relative to WT, suggesting that OsRDRI plays a 
role in genome-wide gene regulation in rice. In addition, 
the Osrdrl mutant showed regional alterations in smRNA 
accumulation and/or titration across the rice genome, and 
at least some of which are associated with locus-specific 
alteration of DNA methylation. 

Among the differentially accumulated smRNAs, many 
were miRNAs both previously known and newly identi- 
fied in this study. Expression changes in many of these 
miRNAs are associated with changes in their target gene 



Wang et al. BMC Plant Biology 2014, 14:1 77 
http://www.biomedcentral.com/1471-2229/14/177 



Page 9 of 12 



Seedling 



Heading 



Maturity stage 




WT 
■ Osrdrl 



CK Sail Cu He NO 



CK S*U Qi Hi NO 




km 



Osrdrl 



Figure 6 Effects of null mutation of OsRDRl on phenotypes under normal and abiotic stress conditions in rice. Phenotyping of Osrdrl 
and WT plants at three growth stages under normal condition and four abiotic stress conditions (salt, heavy metal Cu 2+ , heavy metal Hg 2+ and 
overdose NO). 



expression, which could partly be responsible for the 
gene expression changes observed in the Osrdrl mutant 
relative to WT. The mechanism by which mutation of 
OsRDRl caused changes in miRNA expression in rice 
was not clear. It was found in Arabidopsis that none of 
the RDRs has a direct role in miRNA biogenesis [2,38]. 
However, RDRs could impact miRNA accumulation indir- 
ectly by either affecting miRNA precursor gene expression 
through the TGS or PTGS pathways [1,2,38] or generating 
dsRNAs that compete for DCL1 function that is required 
for miRNA biogenesis [1,2]. 

Our results suggested that OsRDRl might play a role 
in maintaining the intrinsic locus-specific DNA methyla- 
tion patterns, as its mutation caused alteration of methy- 
lation at some of the loci we analyzed. In particular, the 
changes in CHH methylation, which are indicative of de 
novo methylation by the RdDM pathway, showed correl- 
ation with changes in smRNA accumulation. In the re- 
spect of reduced CHH methylation and concomitant 
reduced smRNA accumulation, OsRDRl may functionally 
resemble the Arabidopsis RDR1 and plays a role in the 
21 nt siRNA-dependent non-canonical RdDM pathway 
[20,21]. However, some of the analyzed loci showed in- 
creased CHH methylation that is associated with increased 
smRNA accumulation in Osrdrl. We should caution that 
because we analyzed only 10 loci, the results observed 
may not be extrapolated to global scale. In Arabidopsis, it 
was documented that mutation of RDR1 resulted in near 
complete loss of methylated cytosine of all three sequence 
contexts (CG, CHG and CHH) within the 4,949 CHH hy- 
pomethylation DMRs (differentially methylated regions) 



between drmll2 and WT [20]. Therefore, genome-wide 
methylation analysis (methylome) of the Osrdrl mutant 
will be required to confirm whether OsRDRl plays a simi- 
lar role globally in rice. 

Previous studies in Arabidopsis and Nicotiana have 
defined an established role of RDR1 in plant virus re- 
sponses [24,26,27]. We showed here that Osrdrl exhibited 
no phenotypic differences from its sibling WT plants 
under normal growing condition, but displayed ephem- 
eral phenotypic fluctuations contingent with presence 
of several abiotic stress conditions. This observation, to- 
gether with the enriched GO categories including those 
involved in metabolic process of the differentially expressed 
genes between Osrdrl and WT (Figure 2b), suggest that 
the effects of altered gene expression due to OsRDRl mu- 
tation has been largely canalized under normal condition 
but can be released by certain abiotic stress conditions 
[39], an issue that merits further investigations. Regardless, 
our results suggest that, apart from its established role in 
the production and amplification of exogenous, virus- 
derived siRNAs (vsiRNAs) in infected plants [26,27], the 
rice RDR1 homolog (OsRDRl) might also play a role in 
certain abiotic stress responses, which however may not 
involve stable epigenetic changes in this respect. In this re- 
gard, it should be emphasized again that the genome-wide 
analyses of both smRNA profiles and gene expression in 
Osrdrl were conducted on plants grown under normal 
conditions. Therefore, further studies are needed to con- 
duct the analyses in plants of the mutant and WT under 
both short- and long-term stress conditions. It would also 
be interesting to analyze the progeny of stress-treated 
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Osrdrl plants to investigate if OsRDRl is involved in 
transgenerational inheritance of stress-induced epigenetic 
changes, if they occurred. 

Conclusions 

How RDR1 affects global gene expression and smRNA 
profiles have not been previously investigated in any 
plant species. By analyzing a null mutant of the rice 
RDR1 gene (OsRDRl), we showed that expression of 
more than 1,000 endogenous genes of diverse gene 
ontology (GO) categories were significantly altered in 
the mutant, indicating a functional role of OsRDRl in 
regulating endogenous gene expression in rice. By 
smRNA deep-sequencing, we found that extensive alter- 
ation in smRNA clusters occurred across each of the 12 
rice chromosomes in the mutant, indicating a role of 
OsRDRl in smRNA biogenesis and/or titration in rice. 
We also found that at least some of the gene expression 
changes are correlated with differences in miRNAs. We 
further showed that changes in smRNAs can be concomi- 
tant with locus-specific alteration of cytosine methylation 
primarily of the CHH contexts, thus linking OsRDRl to 
DNA methylation in rice. Finally, we showed that whereas 
no apparent phenotypic abnormality was associated with 
loss of function of OsRDRl, ephemeral phenotypic fluctu- 
ations could be generated by various short-term abiotic 
stress conditions as a result of OsRDRl mutation, suggest- 
ing a role of OsRDRl in plant abiotic stress response. 

Methods 

Plant materials 

Based on information about the rice retrotransposon 
Tosl7 insertion lines (http://tos.nias.affrc.go.jp/), we ob- 
tained a line (#RDR704) of rice cultivar Hitomebore with 
Tosl 7 being inserted into the second exon of OsRDRl in 
a heterozygous state (accession # H0643). According to 
BlastN search at the NCBI website (http://blast.ncbi.nlm. 
nih.gov/Blast.cgi), we found that rice contains a single 
copy of the insert gene (OsRDRl), The three OsRDRl 
genotypes, WT (RDR1/RDR1), heterozygous (RDRl/rdrl) 
and homozygous mutant (rdrl/rdrl) were identified by 
two pairs of specific primers (Figure la). Specifically, WT 
was identified by a pair of primers anchored within the 
OsRDRl gene but flanking the Tosl 7 insertion site; homo- 
zygous mutant was identified by a pair primers with on 
anchored to the OsRDRl gene and the other one targeting 
to the terminal of Tosl 7; and the heterozygote was iden- 
tified by combinations of both types of primers. The 
heterozygotes of OsRDRl(+/-) were selfed for five succes- 
sive generations, and in each of the first four generation 
(S1-S4) only heterozygous individuals were selected based 
on PCR identification. At the last generation (S5), the 
newly segregated homozygous mutant, i.e., OsRDRl-/ - 
(designated as Osrdrl), and its sibling wild type (WT), i.e., 



OsRDRl+l+, were selected and propagated for an add- 
itional generation to have sufficient seeds for this study. In 
this way, the mutant and WT should be genetically identi- 
cal except for the locus in concern, i.e., OsRDRl. Seeds of 
the two genotypes were thoroughly washed with distilled 
water and then germinated in the dark in Petri dishes con- 
taining distilled water at 28°C. After a 2-day incubation, 
germinated seeds were transferred to a greenhouse at 
26°C under 16 h/8 h light/dark regime for the four kinds 
of abiotic stress treatments: 0.15 mMol/L NaCl (salt), 
0.25 mMol/L CuS0 4 (heavy metal Cu 2+ ), 0.25 mMol/L 
HgCl 2 (heavy metal Hg 2+ ), 1 mMol/L Sodium nitroferri- 
cyanide(III) dehydrate (SNP, for NO stress) in Hoagland 
nutrient solution for 7- day. Mock controls (CK) were 
grown in parallel. Then, all seedling plants were trans- 
planted to normal paddy field. Plants were surveyed at 
appropriate growth and developmental stages, seedling, 
heading and maturity. 

Genomic DNA was isolated from seedlings of Osrdrl 
and WT at the same developmental stage using a modified 
CTAB method. Total RNA was isolated from the same 
seedlings with the Trizol Reagent (Invitrogen) according 
to the manufacturers instructions. The RNA was then 
treated with RNase-free DNase I (Invitrogen) to eliminate 
possible genomic DNA contamination before being re- 
verse transcribed with the Superscript RNase H- Reverse 
Transcriptase (Invitrogen). 

SmRNA library construction and sequencing 

Total RNA was prepared for smRNA sequencing based on 
the Illumina Sample Preparation Protocol. The samples 
were quantified and equalized so that equivalent amounts 
of RNA from Osrdrl and WT were analyzed. In brief, total 
RNA was purified by electrophoretic separation on a 15% 
TBE-urea denaturing PAGE gel and smRNA regions cor- 
responding to the 15-30 nucleotide bands in the marker 
lane were excised and recovered. The 15-30 nt smRNAs 
were 5' and 3' RNA adapter-ligated by T4 RNA ligase and 
at each step length validated and purified by urea PAGE 
gel electrophoretic separation. The adapter-ligated smRNA 
was subsequently transcribed into cDNA by Super- Script 
II reverse transcriptase (Invitrogen) and PCR amplified, 
using primers that anneal to the ends of the adapters. The 
amplified cDNA, too, was purified and recovered. The final 
quality of the library was ensured by validation of the size, 
purity and concentration using an Agilent Technologies 
2100 Bioanalyzer. The two constructed cDNA libraries 
subsequently underwent Solexa/Illuminas proprietary 
flowcell cluster generation and bridge amplification. 

Analysis of smRNA clusters 

SmRNA reads of 18-26 nt in size were counted within 
every sliding 100 bp window along the rice genome. The 
reads were normalized to RPM (reads per million), and 
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comparison was then made between Osrdrl and WT 
plants using the median RPM values, which are denoted 
as X for Osrdrl and Y for WT plants. The fold value was 
calculated by the formula log 2 X-log 2 Y = log 2 (X/Y). 

Affymetrix GeneChip® Rice Genome Array 

The microarray transcriptional profiling was performed by 
the Affymetrix, Inc. in the Gene Company Ltd. (Shanghai, 
China), using procedures described in the GeneChip® 
Expression Analysis Technical Manual. 

Real-time quantitative (q) reverse transcriptase (RT)-PCR 
analysis 

The qRT-PCR experiments were performed using SYBR 
Premix Ex Taq (TOYOBO) according to the manufac- 
turers instruction on a Roche LightCycler480 apparatus 
(Roche Inc.). The primers for amplifying the 18 studied 
genes were designed using the Primer 5 software and 
listed in Additional file 6. Primers for qRT-PCR analysis of 
transposase genes were described in a previous report 
[40]. Expression of a rice fi-actin gene, eEF gene and 
UBQ5 gene were used as internal control with the primer 
pairs ofS'-atgccattctccgtctt and S'-gctcctgctcgtagtc; 5'-tttc 
actcttggtgtgaagcagat and S'-gacttccttcacgatttcatcgtaa; 5'- 
accacttcgaccgccactact and 5'-acgcctaagcctgctggtt, respect- 
ively. Conditions of RT-qPCR were as reported [40]. 

Hemi-nested RT-PCR for detecting miRNA expression 

Four miRNA-specific oligonucleotide primers were de- 
signed according to a previous report [41,42] to reverse 
transcribe the specific miRNAs, and four pairs of specific 
primers corresponding to the RT primers were designed 
to amplify the cDNA (Additional file 6). RT-PCR was 
performed as previously described [41,42]. 

Bisulphite sequencing 

DNA (~1 ug) from each plant was treated with bisulphite 
using the EZ DNA Methylation-Gold Kit and amplified 
using specific primers (Additional file 6). Bisulphite PCR 
product was cloned into vector NTI and positive clones 
were sent for sequencing. The results were analyzed using 
analysis software on the website http://katahdin.mssm. 
edu/kismeth based on a previous report [43]. 

Additional files 



Additional file 3: Figure S2. Putative novel miRNAs identified from 
Osrdrl and wide-type, (a) The sequences, expression status including read 
counts, and genomic locations of the novel miRNAs. 5p, the mature 
miRNA sequence resides in 5' half of the predicted stem-loop structure; 
3p, the mature miRNA sequence resides in 3' half of the predicted 
stem-loop structure, mfe, minimum free energy, (b) The predicted 
stem-loop structure of precursor RNA of the novel miRNAs. The mature 
miRNA sequence inside the stem-loop is indicated by a red line, and the 
5' to 3' direction of a miRNA is indicated by an arrowhead. 

Additional file 4: Figure S3, (a) Correlation of expression of known 
miRNAs and their targets, (b) Pairwise comparison between expression 
levels of known miRNAs and their targets. Yellow and blue columns 
represent target and miRNA expression levels, respectively, /-axis 
indicates the values of log 2 fold change. 

Additional file 5: Figure S4. Regional association between smRNA 
clusters and DNA methylation for each of the 10 assayed genomic loci 
from two transposable elements (TEs), TosU (four regions) and Pong 
(three regions) and three genes (one region each). The red, blue and 
green circles denote for CG, CHG and CHH sequence contexts, 
respectively, wherein the filled ones are methylated and empty ones are 
unmethylated. 

Additional file 6: Table S2. Primers for q RT-PCR assay of 18 genes, 
bisulfite sequencing of 10 loci and semi-nested qRT-PCR analysis of four 
miRNAs. 
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